Nuclear "pasta" structures in low-density nuclear matter and neutron star crust 
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Abstract 

In neutron star crust, non-uniform structure of nuclear matter is expected, which is called the "pasta" structure. From the 
recent studies of giant flares in magnetars, these structures might be related to some observables and physical quantities of 
the neutron star crust. To investigate the above quantities, we numerically explore the pasta structures with a fully three- 
dimensional geometry and study the properties of low-density nuclear matter, based on the relativistic mean-field model and 
the Thomas-Fermi approximation. We observe typical pasta structures for fixed proton number-fraction and two of them for 
cold catalyzed matter. We also discuss the crystalline configuration of "pasta". 
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I. INTRODUCTION 

A neutron star has a radius of about 10 km and the 
mass about 1.4 times the solar mass. It is considered 
to consist of four parts [l|, Q- The region around 0.3km 
from the surface is called "outer crust", where Fe nu- 
clei are expected to form a Coulomb lattice. The region 
around 0.3-1 km is called "inner crust" with a density 
about 0.3-0.5 po- There are neutron-rich nuclei in lattice 
and dripped neutrons in a superfluid state. Two cen- 
tral regions with higher densities are called "outer core" 
and "inner core" , where speculated are proton supercon- 
ductivity, neutron superfluidity [3| , meson condensations 
[3-Ql J hyperon mixture @t111 , or quark matter [l^, [l3 ■ 
The transition of matter composition with the change 
of density inside neutron stars causes a question: does 
it change smoothly or suddenly? A sudden change of 
matter property is generally accompanied by a first-order 
phase transition which causes an appearance of the mixed 
phase. Ravenhall et aL[lJ| suggested the existence of 
non-uniform structures of nuclear matter, i.e. the struc- 
tured mixed phase. They suggested five types of geo- 
metrical structures such as droplet, rod, slab, tube, and 
bubble 47] . Due to its geometrical shapes which depend 
on the density, we call it "nuclear pasta" like spaghetti 
and lasagna etc [ij, [l3|- Many studies have suggested 
the existence of the pasta structures in low-density nu- 
clear matter, relevant to the crust region of neutron stars 
and the collapsing stage of supernovae. The existence of 
the pasta structures in the crust of neutron stars may 
not have influence on the bulk property and structure of 
neutron stars so much. However, it should be important 
for the mechanism of glitch, the cooling process of neu- 
tron stars, and the thermal and mechanical properties of 
supernova matter [20| . 

Recently, in the X-ray afterglow of giant flares quasi- 
periodic oscillations (QPOs) have been observed in some 
soft-gamma ray repeaters (SGR) [2l|,|22|- These flares are 
energetic 7-ray bursts from strongly magnetized neutron 



stars, magnetars. As one comprehensible understanding 
of the QPOs they may be attributed to the shear oscil- 
lations of neutron star crust. Based on this interpreta- 
tion, detailed information on neutron star crust can be 
extracted from the QPOs [23|, where the frequencies of 
shear oscillations depend on shear moduli of the Coulomb 
lattice of nuclei. Ogata et al. have calculated shear mod- 
uli using molecular dynamics simulations .24] . However, 
in this calculation, some important effects are missing, 
such as charge screening, finite-size effects of nuclei, su- 
perfluidity of dripped neutrons and so on [25] . Consider- 
ing the realistic situation, these elements should be taken 
into account in the elaborate calculation of shear moduli. 
Pasta structures should be also considered. Thus the en- 
ergy change against a small deformation of the lattice can 
be discussed using shear moduli. On the other hand, if 
the energy change against a large deformation is known, 
the breaking strain of neutron star crust will be obtained. 
Then, this breaking strain could determine the possible 
size of mountains on the neutron star crust, which may 
radiate gravitational waves strongly with rapid rotation 
[23]. These waves may be caught by large-scale interfer- 
ometers and could limit the spin frequencies of accreting 
stars . Furthermore, the breaking strain may be impor- 
tant for the "star quake" model of giant flares in magne- 
tars [131. 

The species and the sizes of the pasta structure are de- 
termined for given average baryon-number densities by 
minimizing the total energy density, which consists of 
the bulk, the surface, and the Coulomb energy densities. 
From the thermodynamical point of view, nuclear matter 
at sub-saturation density can be represented by a dilute 
gas phase and a dense liquid phase in chemical equilib- 
rium which determines particle densities in both phases. 
Once the averaged density is given, the volume fraction of 
gas or liquid phase is determined. Thus, the bulk energy 
density is independent of the shape and the size of the 
structure, which are determined by the balance between 
the Coulomb and surface energy densities [28j . 



In many studies about nuclear pasta, the geometrical 
symmetry of the structure has been assumed by employ- 
ing the Wigner-Seitz (WS) approximation [ij, \^. In 
this approximation, one can find all the physical quan- 
tities from those in a single WS cell. Furthermore, the 
calculation is reduced to the one-dimensional form owing 
to the geometrical symmetry, which drastically saves the 
computational cost. On the other hand, by virtue of the 
recent development in computational science, it becomes 
possible to calculate without any assuruption about the 
geometry of non- uniform structures [29|, 13(| . These stud- 
ies have imposed the periodic boundary condition in a 
small cubic cell, which included only one period of struc- 
tures. Although they got essentially the same results, 
i.e. typical pasta structures, it is hard to extract further 
information such as crystalline configuration of "pasta" 
and mechanical properties. 

For detailed studies about the properties of neutron 
star crust, it is desired to perform three-dimensional 
calculations in a periodic cubic cell with sufficiently 
large sizes. Accordingly, we have developed a numeri- 
cal code to calculate the density distribution of parti- 
cles in a three-dimensional coordinate space and the rel- 
evant physical quantities, e.g., energy density or pres- 
sure of matter, by using a relativistic mean-field (RMF) 
model under the Thomas-Fermi approximation [31]. We 
have explored ground states of low-density nuclear mat- 
ter with a fixed proton number-fraction and observed 
a series of typical "pasta" structures appearing as the 
ground states. We have also observed crystalline con- 
figurations of droplets, rods, tubes, and bubbles. One 
of our findings is the appearance of face-centered cubic 
(fee) lattice of droplets in the ground state. However, 
the system was limited to the cases with fixed proton 
number-fraction and detailed discussion is still needed 
about its appearance. 

In this article, we explore the pasta structures and 
properties of low-density nuclear matter, following the 
line of the previous study in Ref. [3l| but in more de- 
tail and with wider scope to include catalyzed matter. 
Appearance of the pasta structures and other unusual 
structures, crystalline configuration of droplets are dis- 
cussed. 



In Sec. II, we present the model and describe our nu- 
merical procedure. In Sec. Ill, we first demonstrate some 
results from our three-dimensional calculation for low- 
density nuclear matter with fixed proton number-fraction 
Yp = 0.5, 0.3, 0.1 that may be related to the super- 
nova matter and newly born hot protoneutron star crust. 
Then, in the second part of Sec. Ill, we investigate the 
pasta structures at /3-equilibrium, as they occur in cold 
neutron stars, and the crystalline configuration of nu- 
clear pasta. Finally, Sec. IV is devoted for summary and 
concluding remarks. 



II. MODEL AND METHOD 

There have been used the many-body techniques for 
studying the pasta structures in the literature, such as 
compressible liquid-drop modeHCLDM) [H, [11, IH-Sl , 
the Thomas- Fermi modelfH, [sB-ill, the Hartree-Fock 
(HF) approximation employing effective NN interactions 



391 . |40|, and quantum molecular dynamics (QMD) 
|43| . In the studies using CLDM and the Thomas- 
Fermi model, always used was the WS approximation, 
where only the typical pasta structures are assumed. 
QMD calculation does not assume any specific non- 
uniform structure of baryon matter, while uniform back- 
ground of electrons, though their distributions should 
be almost uniform, is assumed. Some of the Thomas- 
Fermi calculations and the HF calculation used a periodic 
boundary condition and did not assume any geometrical 
symmetry for the structure. However, the size of the 
periodic unit cell was not large enough for quantitative 
discussion. In this paper, to describe interaction among 
nucleons, we employ an RMF model under the Thomas- 
Fermi approximation [38|. It deals with the mean- fields 
of sigma, omega, and rho mesons represented by u,uj^^^ 
and R^, respectively, together with nucleons i/', electrons 
■0e, and electromagnetic field A^, by a Lagrangian C in- 
troduced in a Lorentz-invariant form as follows. 
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where U{a) — ■^brnpflg^Na)^ — jclgaNcr)^ is a non-linear 
term for the scalar field, m^(r) = toa? — g^NC^i'^) rep- 
resents an effective mass of nucleon, and gab represent 
coupling constant between a and b. Thermodynamic po- 
tential n is then given hy Vl = E - J2i=„,p,e f^i I d^'^Pii^) 
in terms of total energy E, number densities Pi{v) and 
chemical potentials /i^. Variation of Q. with respect to 
each field gives the field equation. 

The set of coupled field equations for the mean-fields 
and the Coulomb potential A^^ — (Vcoui,0), render 

- V2(7(r) + mla{Y) = g,M (p;(r) + p,^,(r)) - —,(2) 

- V^woCr) + m^ wo(r) = 5^Ar(pp(r) + p„(r)), (3) 

- V2i?o(r)+m2i?o(r)=.gpA,(pp(r)-p„(r)), (4) 
Wcoui(r) = 47re2 (pp(r) - p,{y)) , (5) 

where pl{v) — {■ij;i{r)'tlji{r)),i = p, n is the nucleon scalar 
density and Rq is the third component of the isovector 
field Rp. Field equations for fermions simply yield the 
standard relations between the densities pi and chemical 



potentials /i^ within the Thomas-Fermi approximation, 



(a) 



M« = Y fcF,„(r)2 + m*j^{rf + 5^jvWo(i") - 5pAri?o(i") 

(6) 

fJ-p = ykF^p{r)^+m*j^{rf + gajN^a{r)+gpNRa{r) 

-"ttoui(r) (7) 
Pe(r) = (Aie-V^Coui(r))V37r2, (8) 

where the local Fermi momentum, kp^i{v), is simply re- 
lated to the density, fc|, j(r)/(37r^) — pi{r). Finally the 
baryon-number conservation and the charge neutrality 
are imposed besides these equations. We use the same 
set of parameters as in Ref. [38] listed in Table |TT1 in 
order to compare the equation of state (EOS) and struc- 
tural changes of the pasta structure with and without the 
WS approximation. With these parameters, we can re- 
produce the properties of uniform nuclear matter shown 
in Table HI The first and second quantities, po and eo, 
are the saturation density of symmetric nuclear matter 
[fv 0.16 fm~'^) and its energy per nucleon, respectively. 
The third and forth quantities K and 5*0 are the incom- 
pressibility and symmetry energy at po ; respectively. The 
last one L is the slope parameter of symmetry energy at 
Pq. Using these parameters the binding energy per nu- 
cleon around the saturation density is expressed as 
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PO [fm-3] 


60 [MeV] 


K [MeV] 


So [MeV] 


L [MeV] 


0.153 


-16.4 


240 


33.4 


84 



TABLE I: EOS of uniform nuclear matter. 



To numerically simulate the non-uniform structure of 
infinite matter, we use a cubic cell with a periodic bound- 
ary condition. If the cell size is small that includes only 
one or two units of the structure, the geometrical shape 
should be affected by the boundary condition and the 
appearance of some structures is implicitly suppressed. 
Therefore, the cell size should be so large as to include 
several units of the pasta structures. We divide each cell 
into three-dimensional grids. The desirable grid width 
should be so small as to describe the detailed density 
distribution, particularly at the nuclear surface. Due to 
this requirement, we set the grid width to 0.3 fm at the 
largest. 

Giving the average baryon-number density pB , initial 
density distributions of fermions are randomly prepared 
on each grid point. Then proper density distributions 
and the meson mean-fields are searched for. We intro- 
duce the local chemical potentials /ia(r) (a = p,n,e) to 
obtain the density distributions of baryons and electrons,. 
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FIG. 1: (a) Schematic figure of density distribution and 
the chemical potential subtracting a constant mass. In the 
region with zero density, where the kinetic energy is zero, 
the chemical potential p{v) is identical to the potential U(v), 
while the chemical potential in the region with finite particle 
density is constant. 

(b) Schematic figure of non-uniform matter with an isolated 
matter region, where the chemical potential is different from 
that of the global one. By the procedure to adjust densities 
between only the neighboring grids, such unphysical region of 
isolated matter may appear. 

The equilibrium state is eventually determined so that 
the chemical potentials are independent of the position. 
An exception is the region with no particle density, where 
the chemical potential of that particle can be higher (see 
Fig.[IJa) for more explanation). We repeat the following 
procedures to attain uniformity of the chemical poten- 
tials. A chemical potential pi{v) of a baryon i = p,n on 
a grid point r is compared with those on the six neighbor- 
ing grids r' = r + dr, {dr = ±dx, ±dy, ±dz). If the chem- 
ical potential at the point under consideration is larger 
than that of another /ii(r) > /ii(r'), some part of the 
density will be transferred to the other grid point. This 
adjustment of the density distribution is simultaneously 
done on all the grid points. In addition to the above pro- 
cess, we adjust the particle densities between distant grid 
points chosen randomly in order to eliminate regions with 
different pi which can happen to isolated matter regions 
as Fig. [IJb) . The meson mean-fields and the Coulomb 
potential are obtained by solving Eqs. ©-(El) using the 
baryon density distributions pi (r) (i — p, n) and the 
charge density distribution ^^(r) — Pe(r). The electron 
density pe (r) is directly calculated from the Coulomb po- 
tential Vcoui(r) and the electron chemical potential /ie, 
Pe(r) = {Pe — Vcoui)^/(37r^). The global charge neutral- 
ity is then achieved by adjusting p^- Above processes are 
repeated many times until we get convergence. 

We have used PRIMERGY BX900 of JAEA massively 
parallel computing system. To complete a calculation of 
a typical case, about 510 CPU hours is needed at least. 



III. RESULTS 



A. Fixed proton number-fraction 



First, we present here some results for fixed proton 
number-fraction Yp with 1^=0.5 (symmetric nuclear mat- 
ter), 0.3, and 0.1 which are roughly relevant to the su- 
pernova and the neutron star crust. Shown in Fig. [5] are 
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m^ [MeV] 


m„ [MeV] 


mp [MeV] 


6.3935 


8.7207 


4.2696 


0.008659 


-0.002421 


400 


783 


769 



TABLE II: Parameter set used in RMF model. 



the proton density distributions in cold symmetric mat- 
ter. We can see that the typical pasta phases with rod, 
slab, tube, and bubble, in addition to spherical nuclei 
(droplet), are reproduced by our calculation in which no 
assumption on the structures was used. Furthermore, 
these cells include several units and we can specify these 
lattice structures. Crystalline configuration of droplet 
and bubble is "fee" , rod and tube is "honeycomb" . 




FIG. 2: (color online) Proton density distributions in the 
ground states of symmetric matter {Yp =0.5). Typical pasta 
phases are observed: (a) Spherical droplets with a fee crys- 
talline configuration at baryon density pB = 0.01 fm""'. (b) 
Cylindrical rods with a honeycomb crystalline configuration 
at 0.024 fm-^ (c) Slabs at 0.05 fm'^ (d) Cylindrical tubes 
with a honeycomb crystalline configuration at 0.08 fm~^. (e) 
Spherical bubbles with a fee crystalline configuration at 0.09 
fm-^ 



Any exotic mixture does not appear as a ground state 
at any density. In a droplet, we have seen that the proton 
density is highest near the surface due to the Coulomb 
repulsion, while the neutron density distribution is flat 
inside the droplet. Note that baryon density outside the 
droplets is zero for 1^=0.3 and 0.5. Electron density is 
spread over all space but slightly localized around the 
droplets which brings about the charge screening effect. 

We can see the density distribution of fermions for 
Yp = 0.5 and ps = 0.03 fm^^ in Fig. |3] and find 
the unique features of the three-dimensional calculation. 
Here, in the upper panel, we show the proton density 
distribution on a sliced plane and depict the two kinds 
of the density profiles of proton, neutron, electron and 
the Coulomb potential in the lower panel. One is along a 
vertical line (a) which passes through the rods, another 
is along a horizontal one (b). From this figure, the ad- 
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FIG. 3: (color online) Upper: Proton density distribution in 

the rod phase {pB ~ 0.03 fm~^] for Yp =0.5) on the sliced 

plane. 

Lower: the density distribution along a line (a) (vertical) and 

(b) (horizontal) in upper figure. Red lines indicate proton, 

green neutron, blue electron and pink the Coulomb potential. 

vantage of the three-dimensional calculation can be seen. 
The proton, neutron, and electron density distributions 
are almost the same for both the cases of (a) and (b). 
However, a slight difference appears in the Coulomb po- 
tential. We have set the maximum value of the Coulomb 
potential to be zero for convenience. In the case of the 
rod phase, that point corresponds to the centroid of the 
triangular lattice. These points are included not on the 
path (a) but on the path (b). Considering the impor- 
tance of the distinct relation between the Coulomb and 
surface energies for the pasta structure, we should take 
into account this anisotropy in a proper way. 

In Fig. |4] we show the energy, total pressure, and 
baryon partial pressure as functions of density. Baryon 
partial pressure is given by subtracting the electron con- 
tribution from the total pressure. Note that the energy 
E/A — rriN includes the kinetic energy of electrons, which 
makes the total pressure positive. This density depen- 
dence is qualitatively the same as the one with the WS 
approximation. The appearance of non-uniform struc- 
tures will make nuclear matter more stable: the energy 
per baryon gets lower up to about 15 MeV/ A compared 
to uniform matter and the pressure per baryon higher up 
to about 0.5 MeV/ A. 



We have obtained almost the same EOS with that 
given by the WS approximation, in which the same RMF 
model is applied. However, one of the differences between 
our results and those with the WS approximation appears 
in the existing region of each pasta structure; the density 
region of the rod is wider and the tube narrower in our 
calculation. Since the energy differences between differ- 
ent structures are quite small, the crystalline configura- 
tion might affect the appearance of each pasta structure. 

Figure [5] shows the radius of droplets Rd, the lattice 
constant i?iattj ^md the volume fraction u in the droplet 
phase. Here, Rd and _Riatt are defined as follows, 
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(11) 



and u — {Rd/RiattY, where V denotes the cell volume, 
Nd the number of droplets in the cell, and the bracket (...) 
means the average over the cell volume. Comparing the 
present results with those using the WS approximation, 
the volume fraction exhibits the same behavior, but the 
lattice constant and the radius of the droplet the differ- 
ent behavior. This might be caused by some differences 
of the treatment of the Coulomb interaction between the 
present calculation and the one using the WS approxi- 
mation. 

Although any unusual structure does not appear in 
the ground state, one new feature emerges in the crys- 
talline configuration of droplets. It has been considered 
within CLDM that a body-centered cubic (bcc) lattice 
is energetically favored than the face-centered cubic (fee) 
lattice by the Coulomb energy [321 • However, in our cal- 
culation, the fee lattice emerges as a ground state at some 
density region. Crystalline configurations of the bcc and 
fee lattices give rise to a subtle difference of the Coulomb 
energy, which amounts to about 0.2-0.8 MeV. The ratio 
of the Coulomb energy to the total energy difference is 
about 20%. In other words, the major contribution to the 
energy comes from the mean-fields and the kinetic ener- 
gies of nucleons, i.e., the bulk energy of nuclear matter. 
Therefore, the ansatz used in the previous studies, e.g., 
CLDM, should be dangerous since small differences in 
size and composition of droplets give a significant change 
in the total energy of the system. Table IIIII shows the 



PB [ fm--'] 


0.012 


0.014 


0.016 


0.018 


0.020 


Rd{icc) [fm] 
7?d(bcc) [fm] 


6.86 
6.99 


7.04 
7.18 


7.23 
7.36 


7.61 

7.75 


7.79 
7.92 



TABLE III: Density dependence of the droplet radii for the 
fee and bcc lattices. 

density dependence of radii of the droplets in the cases 
of the fee and bcc lattices. The radius of the droplet is 
different even if their baryon densities are the same. In 
Refs. ^j, it has been reported that the bcc lattice is 



realized in the ground state, by comparing the bcc and 
fee crystals of droplets of the same size, while they have 
different droplet sizes in our case. 

In the QMD calculations [42] that precede the present 
calculation without assuming geometrical structures, 
droplets form the bcc lattice. This difference might come 
from the treatment of electrons or the charge screening 
effect: uniform electron distribution has been assumed 
in the QMD calculation. To see the effects of the elec- 
tron distribution on the crystalline configuration, let us 
compare two cases within our framework: one is the full 
calculation and the other assumes uniformly distributed 
electrons. In the latter calculation, the Coulomb poten- 
tial Vcoui in Eq.® is replaced by a constant Vb = and 
Pe = {fie — Vb)'^/37r^. In the full calculation using ([5]), 
the value of Vq is arbitrary due to the gauge invariance, 
and one can take Vq for the sake of convenience, e.g. , 
as Vo =0. However, in the case of uniformly distributed 
electrons, the gauge invariance is partially violated, since 
we replace Vcoui to Vq in the equation for the electron 
chemical potential but retain Vcoui in the equation for 
the proton chemical potential and thus in the expression 
for the proton number density. In this case, the droplets 
have a smaller size compared with that in the full cal- 
culation in all the density region of the droplet phase. 
This means that the Coulomb repulsion among protons 
is slightly weaker in the full calculation due to the screen- 
ing by electrons. One can see this difference in Table 
IIVI In our previous study using the WS approximation 



Pb [fm 



0.004 



0.01 



0.016 



0.022 



Rd (full calc) [fm] 
Rd (unif clec) [fm] 



6.09 


6.67 


7.23 


7.79 


5.91 


6.40 


7.01 


7.64 



TABLE IV: The radii of droplets by the full calculation and 
with uniformly distributed electrons. 



[38J, a similar discussion was made. We argued that the 
charge screening effect is not so remarkable due to the 
large Debye screening length. However, the crystalline 
configuration does not change and the fee lattice is more 
favored. Thus, it is confirmed that the charge screen- 
ing by electrons does not very much affect the crystalline 
configuration. This difference of the crystalline config- 
uration between the QMD calculation and the present 
calculation remains to be elucidated as a future problem. 
We obtain the typical pasta structure as a ground state 
for any proton number- fraction above 0.1. Also in the 
cases of 1^=0.3 and 0.1, the fee lattice of droplets is en- 
ergetically more favorable than the bcc one. In Fig. [6] 
we depict the density profiles of proton, neutron, elec- 
tron and the Coulomb potential for Yp — 0.3 and 0.1 
with baryon number density 0.03 fm~'^ along a line which 
passes through the rods in the same way as in Fig. [S] 
While for Yp =0.3 there appear vacant regions of neutron, 
the neutron density is finite at any point for Yp =0.1: the 
space is filled with dripped neutrons. Even in the case of 
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FIG. 4: (color online) From the left, energy, pressure, and baryon partial pressure of Yp=Q.^, 0.3, 0.1 in this order from the 
upper panel. Red lines indicate droplet, green rod, blue slab, magenta tube, cyan bubble, and black uniform, respectively. 
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FIG. 5; (color online) Density dependence of the radius, the 
lattice constant, and the volume fraction. Red lines are our 
results and green lines are those with the WS approximation. 

1^ = 0.1, we can see that the proton density is highest 
around the surface due to the Coulomb repulsion. 



B. Catalyzed matter 

Cold catalyzed matter requires beta equilibrium in- 
stead of the fixed proton number- fraction: fin = fJ-p + fJ-e- 
Shown in Fig. [7] are the proton density distributions of 
the ground states in cold catalyzed nuclear matter. We 
have obtained the body-centered cubic (bcc) lattice of 
droplets, face-centered cubic (fee) lattice of droplets, and 
honeycomb lattice of rods, depending on density. In 
FiglSl the density profiles of fermions in the bcc and fee 
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FIG. 6: (color online) Density distribution of proton, neu- 
tron, electron and the Coulomb potential of rod in ps = 0.03 
fm~'^ (upper; Yp — 0.3, lower: V^ = 0.1) Red lines indicate 
proton, green neutron, blue electron, and pink the Coulomb 
potential. 



lattices of droplets and the Coulomb potential are de- 
picted along a line which passes through the droplets for 
PB ~ 0.01 fni-3 and ps « 0.03 fm-^. The effect of the 
Coulomb repulsion can be seen where the proton density 
is highest near the surface, while the neutron density dis- 
tribution is flat in the droplet. We have observed only 
those three structures (a), (b), (c) in Fig. [T] This result 
is consistent with the previous study about the relation 
between the density region of the pasta structure and 
the slope parameter L in Eq. (10) 44]. The larger L is 
given, the narrower pasta region is. In our calculation, L 
is about 80 MeV, which is close to the critical value « 90 
MeV in the CLDM, where the pasta structures do not 
appear. In the previous calculation using the WS approx- 
imation, only the droplet structure appeared as a ground 
state |38j. However, in our three-dimensional calculation 
using the same RMF framework, the rod structure also 
appears. 
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FIG. 7: (color online) Proton density distributions are shown; 
droplet (bcc) at ps = 0.01 fm" , droplet (fee) at ps ~ 0.03 
fm""^, rod (honeycomb) at ps — 0.056 fm""^. 
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FIG. 8: (color online) Density distributions of proton, neu- 
tron, and electron and the Coulomb potential in the fee and 
bcc lattices. The upper panel shows the case of ps — 0.01 fm^^ 
and pB = 0.03 fm~^ in the lower panel. Red lines indicate the 
proton, green neutron, blue electron (50 times) density distri- 
butions, pink the Coulomb potential, and black dashed lines 
are the lattice constant. 



We show the density dependence of the total energy, 
the Coulomb energy and the proton number-fraction in 
Fig. [9l To see the difference between droplets in the bcc 
and fee lattices, we plot the density dependence of the 
size of droplet, lattice constant, volume fraction, proton 
number of each droplet and proton number-fraction in 
Fig. \T0\ Here, the meaning of Rd and i?iatt are the same 
as in Eqs. dTU]) and pl|) . 

The density dependence of the ground state energy is 
shown in the left panel of Fig. [9l We can see that the 
ground state configuration changes depending on density. 
In the lower density region, the bcc lattice of droplets ap- 
pear. Around pB ~ 0.01 fm~^, lattice structure changes 
from the bcc to fee lattice. Remarkable change occurs 
around pb « 0.052 fm~^: from the fee lattice of droplets 
to the honeycomb lattice of rods. It is hard to see 
the total energy difference between the bcc and fee lat- 
tices of droplets. That of the proton number-fraction at 
PB < 0.01 fm^'^ is also hard to be distinguished in the 
right panel of Fig.[9l However, there are significant differ- 
ences in the proton number-fraction at pB > 0.01 fm~'^ 
and the Coulomb energy as shown in the middle panel of 
FigEl 

As in the case of the fixed proton number-fraction, 
there emerges a fee lattice of droplets near the transition 
density from the droplet phase to another in our calcu- 
lation, while it has been regarded to take a bcc lattice 
in the previous studies [32| . Almost the same radii and 
density distribution of droplet, proton number- fraction 
and proton number in nuclei are obtained for both crys- 
talline configurations (see the right panel of Fig. [TU] at 
Pb < 0.01 fm"'^). The difference between the bcc and 
fee lattices may be seen only in the Coulomb energy: 
the Coulomb energy of the fee lattice is a little higher 
than that of the bcc lattice. However, near the transi- 
tion density from droplet to rod, the radius of droplets 
and the proton number-fraction are different between the 
bcc and fee lattices even if their baryon-number densities 
are the same (see Fig. [TU]at pB ~ 0.02 fm~'^): the size 
of the droplet and the proton number-fraction in the fee 
lattice is 7.54 fm and 0.016, respectively, while those in 
the bcc lattice are 7.01 fm and 0.014. Because the size 
of droplets and the proton number-fraction are different, 
the Coulomb energy alone is no more the criterion of the 
ground state. We should take into account the size of 
droplets and the proton number-fraction looking for the 
ground state in a self-consistent way. Roughly speaking, 
the larger the radius of droplets is, the smaller the surface 
energy is. While the Coulomb energy of the fee lattice 
is larger than that of the bcc lattice in all the region of 
the droplet phase, the total energy of the fee lattice is 
less than that of the bcc lattice by the gain of the surface 
energy. 

Because we cannot take into account the shell effects 
like Ref. [431, the proton number continuously decreases 
with increase of baryon number density. To compare 
the density dependence of the proton number within the 
same interaction for baryons, our model might be simi- 
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FIG. 9: (color online) From the left, the total energy, the Coulomb energy, and the proton number-fraction. Red lines indicate 
droplet (fee), green droplet (bcc), blue rod, respectively. Dotted lines are for the case of a single phase. 
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FIG. 10: (color online) In the left panel, the density de- 
pendence of radii of nuclei, lattice constants, and the volume 
fraction are shown. The proton number of nuclei and the pro- 
ton number-fraction are in the right panel. The lines in each 
panel indicates fee by the red line, bcc by the green line in 
comparison with that in the WS approximation. 



by using a uniform background electrons and uniform 
baryon density inside a nucleus. The effect of the screen- 
ing by charged particles, which is properly included in 
our calculation, may be one of the origins of this differ- 
ence. To get the final conclusion of this issue, we should 
perform another calculation with uniformly distributed 
electrons and confirm the effects of non-uniformity of 
electron. Also the difference of the droplet surface may 
give another reason: they used the CLDM with a sharp 
surface, while our droplets have a diffuse surface. The 
same discussion may apply to the QMD calculation for 
the case of the fixed proton number-fraction. 

At low density, where the volume fraction of droplet is 
less than 1%, we can regard nuclei as point particles even 
if there is a broad distribution of neutron. However, near 
the transition density from droplet to rod, it amounts to 
about 10% and the distance between the nuclei is very 
close to each other, which is one reason of the shape 
transition. In this region, the approximation for each 
droplet as a point particle is not allowed. 



lar to the type B in Ref. [4J|. In our model, the slope 
parameter of symmetry energy L is a little larger than 
that in Ref. [4J|, so that, in our model, the proton num- 
ber converges to about 35 in the low-density limit, while 
about 40 in Ref. [4J|. This result is consistent with the 
relationship between L and the proton number in the 
droplet. 



The maximum size of the droplet can be estimated by 
the Bohr- Wheeler condition as £;*?j^j > 2Es^,f ^. The 
"virial" theorem for the pasta structure reads i^surf = 

2£'coui where i^coui i^ ^^"^ Coulomb energy of an isolated 
nucleus and the Coulomb energy of a nucleus in matter 
as ii^coui ~ ^Coui(l ~ 3u^/'^/2). From these equations, 
the appearance of non-spherical nucleus in nuclear mat- 
ter has been expected for the volume fraction u > 0.125. 
However, in our calculation, the structural change from 
droplet to rod occurs around u « 0.1 (see the value at 
Pb ~ 0.05 fm~'^ in Fig.fTO]). This means that considering 
the non-uniformity of electron is worthwhile for "pasta" 
structures, because the relation between the Coulomb en- 
ergy of a cell and that of a nucleus has been derived 



IV. SUMMARY AND CONCLUDING REMARKS 

We have numerically explored non-uniform structures 
and discussed the properties of low-density nuclear mat- 
ter and cold catalyzed matter, using the RMF model 
under the Thomas-Fermi approximation. Without any 
assumption about the geometric structure, we have car- 
ried out fully three-dimensional calculations in large cu- 
bic cells with the periodic boundary condition. 

First, we explored the low-density nuclear matter with 
the fixed proton number-fraction Yp = 0.5, 0.3, 0.1, 
which may be relevant to supernova explosions and newly 
born protoneutron stars. With increase of density, which 
ranges from well below to around the normal nuclear 
density, we have observed the typical pasta structures 
as a ground state for each density and proton number- 
fraction. The appearance of the pasta structures lowers 
the energy, while the energy differences between various 
geometrical structures are very small. More exotic struc- 
tures like "diamond" , "dumbbell" , and mixtures of two 
types of pasta structures appear as metastable states at 
some transient densities. 



Secondly, we extend to the cold catalyzed matter which 
corresponds to the neutron star crust. In this case, with 
increase of density, which ranges from well below to half 
of the normal nuclear density, we have observed that the 
ground state of matter shows two types of pasta struc- 
tures, droplets and rods. As for the crystalline configu- 
ration of droplets, near the transition density to rod, the 
fee lattice is more favorable than the bee lattice, which 
is different from the previous studies. We have discussed 
some reasons of the difference, but more elaborate studies 
are needed to clarify it. 

In this article, the ground state of low-density nuclear 
matter changes the form of crystalline configuration from 
the bcc to fee lattice near the half of normal nuclear 
density within the RMF model. This conclusion arises 
from including the smooth surface of nuclei and the self- 
consistent calculation of the Coulomb interaction. On the 
other hand, there are several forms of lattice structures 
like hexagonal closest packing, tetrahedron. The hexag- 
onal lattice has the same filling factor with the fee one. 
There are many possibilities, but is not specific inter- 
pretation for the crystalline configuration in the ground 
state. 

In application to neutron star crust, it is interesting 
and might be important to investigate the shear mod- 



ulus. Using the density distribution obtained by our 
three-dimensional calculation for cold catalyzed matter 
(droplet and rod), we can estimate the shear modulus, 
in which the charge screening effect and the finite-size 
effects are properly taken into account. More interest- 
ingly, one may directly calculate the breaking strain by 
considering the large deformation of the lattice. These 
issues will be discussed in a separate paper. 

For newly born neutron stars, as in supernova explo- 
sions, finite temperature and neutrino-trapping effects 
become important, as well as the dynamics of the first 
order phase transition with formation of the structures. 
It would be also interesting to extend our framework to 
include these effects. 
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